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Abstract 

We propose a new feature extraction method based on two dynamical systems induced by intensity 
landscape: the negative gradient system and the Hamiltonian system. We build features based on the 
Hamiltonian streamlines. These features contain nice global topological information about the intensity 
landscape, and can be used for object detection. We show that for training images of same size, our 
feature space is much smaller than that generated by Haar-like features. The training time is extremely 
short, and detection speed and accuracy is similar to Haar-like feature based classifiers. 

1 Introduction 

Feature extraction is a critical component in pattern recognition and image processing. In image pattern 
recognition, one usually designs a family of features, which provides a mapping from image space to feature 
space. If the features are well designed, one can expect that images from different classes of objects (say, 
human faces and cars) are well separated in the feature space. Abundant types of features exist, e.g. Haar- 
like features |[T7l , interest points operators [9 ], template matching |fl5l , etc. However, many of these features 
are computationally expensive. For example, one popular method for object detection is design a large 
family of Haar-like features and use boosting method for feature selection. Consider an n x n training 
image, the number of all possible Haar-like features is of order f2(ra 4 ). In practice, one only affords to 
select a fraction of these features, so the resulting subset is still of order Q(n 4 ) asymptotically. Since the 
complexity of a boosting method is proportional to the number of features, training process is usually slow, 
despite the fact that we can evaluate single feature value fast by using integral image. 

In this paper, we propose a new feature extraction methodology. Our philosophy is that in many cases, 
we can do object detection by inspecting shape characteristics of the intensity landscape. For example, in 
a typical image landscape of human face, eyes and mouth can be identified as valleys, nose and cheeks as 
peaks. If we can find an efficient way encoding this landscape, we can, up to diffeomorphisms, recover 
the whole intensity landscape by decoding. Y.Shinagawa [14] proposed a surface coding method based on 
Morse theory, which represents a closed 2-dimensional manifold by its Reeb graph. A major result in Morse 
theory states that the global shape of a closed manifold is completed determined (up to diffeomorphisms) 
by local critical points of a Morse function. This provides a powerful tool to infer global information from a 
few local information. However, Morse theory only deals with topological invariants, thus, much geometric 
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information is lost in the resulting Reeb graph. (We refer readers to [10] for a complete introduction on 
Morse theory). In this paper, we propose a method that can capture both topological and geometric features 
of such landscapes. Our method is based on two dynamical systems induced by the intensity landscape: a 
negative gradient flow system and a Hamiltonian system. By applying results in dynamical systems, we are 
able to capture the characteristics of the surface by a family of contours and the gradient information along 
the contours. It turns out that, empirically, there are about 0(n) such features for n x n training images. Our 
experiments on face detection shows that one can build an AdaBoost classifier ifTTl based on these features 
and the resulting classifier has high testing accuracy. Most importantly, since the feature space is much 
smaller comparing to that of Haar-like features, training time is significantly shorter. 

It is worthwhile to note that our feature extraction method, when applied to face detection, is far different 
from existing methods. Several popular existing methods are: 1) template matching lfT8l l6l: 2) PCA based 
methods lfl6ll : 3) Neural network based methods [ 12]; 4) Haar-like features with AdaBoost algorithm [17]. 
See [13 ] for a survey. The first type of method typically requires manually designed templates or large 
amount of manually labeled landmark points. The last three types of methods, although fully automated, 
are computationally intensive. Our algorithm has merits from both sides: it can intelligently capture a few 
face related feature candidates and do feature selection on a relatively small feature space, therefore achieve 
very fast automatic training. 

The rest of this paper is organized as follows: in Section 2 we review some background of dynam- 
ical system theory. In Section 3 and 4, we describe algorithms for feature extraction. This is the main 
contribution of this paper. In Section 5 we show the general framework of our training process and some 
experimental results. In Section 6 we discuss limitations of our method and possible future work. 

2 Review of Dynamical Systems 

Any smooth function defined on the plane can be visualized as a landscape. This landscape naturally induces 
a gradient flow system, whose solution curves reflect geometric properties of the underlying landscape. 
An image can be viewed as a discretized landscape, which admits a discretized gradient flow. Our basic 
philosophy is to actively capture features in the image by studying geometric and topological properties of 
the flow. First, we provide a brief review of related topics in dynamical systems. For more information, we 
refer readers to [5 ] and H. 

2.1 Gradient System and Hamiltonian System. 

A planar dynamical system is a system of differential equations: 



where x € R , y € R, f(x, y) and g(x, y) are continuous real-valued functions. The dot means differentiate 
with respect to time t. A curve -f(t) = (X(t),Y(t)) is a solution curve of ©, if X(t) = f(X(t),Y(t)) 
and Y(t) = g(X(t),Y(t)) for any t where j(t) is defined. Note (Q]) defines a velocity field on plane, which 
generates a flow : R 2 — > R 2 , where 4>t(xo, yo) can be obtained by starting at point (xo, yo) and follow the 
velocity field for time t. If we replace f(x, y) and g(x, y) by —f(x, y) and —g(x, y) respectively, we get a 
backward flow of (Q]). (Q]) is called a gradient system if there is a function H (x, y) such that f(x, y) = ^ 
and g(x, y) = ^M- . Conversely, given H (x, y), one can define a negative gradient system 
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(c) Negative Gradient Flow (d) Hamiltonian Flow 

Figure 1 : Example of intensity landscape and the two dynamic systems 



Dynamical behavior of Q is straightforward: at each point, the vector field points to the direction where 
H[x, y) decays fastest. Therefore, local minima of H become sinks in the system, whereas local maxima 
become sources. 

Given a gradient flow system, there exists a Hamiltonian system associated with it. For the same 
H(x, y), one can define 

x 



dH 
dy 

y dx 



(3) 



An important property is that solution curves of ([3]> are exactly level curves of H. To see this, let (x(t),y(t)) 
be any solution, we have 

d ,„, , s , ... dH dx dH dy 
-( ffWt)>!/ ( ( ),, = _- + _-|.0 V( 

It is also worthwhile to note that solution curves of the gradient system and those of the Hamiltonian system 
form an orthogonal family, because the two vector fields are perpendicular to each other at every nonsingular 
point. 

In Figure [U we show an average face image, its intensity landscape, and the associated negative gradient 
flow and the Hamiltonian flow. 



2.2 Singular Point, Attractor and Repellor 

A set is an invariant set if it is a union of solution curves. Singular points and periodic orbits are important 
invariant sets in a dynamical system. Typical singular points in a gradient system are sinks, sources, and 
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saddles. Typical singular points in a Hamiltonian system are centers and saddles. Attractors and repellors 
are generalized concepts of sinks and sources. A region R is attracting if for any point x G dR, x ■ t ^ R 
for any t < 0. (Here, x ■ t means start with point x, follow the flow for time t). Similarly, R is repelling if 
for any point x G dR, x ■ t £ R for any t > . Intuitively, for intensity landscape, one can identify valleys 
as attractors and ridges as repellors. In fact, ridge detection has been discussed in literature extensively. See 
iTSTl for example. However, one limitation of direct ridge detection is that it is highly sensitive with respect 
to random perturbation. In 00, scale-space theory is used to overcome some problems. 

In our method, we do not explicitly look for the exact locations of these invariant sets. Instead, we look 
for certain curves or neighborhoods that contain these invariant sets. By using tools in dynamical systems, 
one can predict the existence and the type of invariant sets inside the curves or neighborhoods. Poincare 
index and Conley index are two nature choices for this purpose. These indexes are invariant under certain 
deformation and perturbation of the underlying dynamical systems. We shall review basic results of these 
two indexes below. 

2.3 Poincare Index 
Definition 2.1. Let 

f x = f(x,y) 

1 y =g(x,y) 

be a planar system. Let C be a simple closed counterclockwise curve that not passing through any singular 
point. Then the Poincare Index of C is defined as 

— / d{arctan(-^-)\ 

The Poincare index has a simple geometric interpretation: consider the orientation of the vector field at 
a point p € C. Let p traverse C, the vector (f(x, y),g(x, y)) rotates continuously and, upon returning to 
the original position, must have rotated through an angle 2irk for some integer k, and this k is the Poincare 
index. Cf. Figure [2] 

The following theorem shows how Poincare index can be used to detect and classify singular points 
enclosed by a curve. See HI for a proof. 



Theorem 2.2. 

(i) The index of a sink, a source or a center is +1 

(ii) The index of a hyperbolic saddle point is —1 
(Hi) The index of a closed orbit is +1 

( iv ) The index of a 

(v) The index of a 



In figure\2\ the vector field is deformed. However, 



closed curve not containing any fixed points is 0. 

close curve is equal to the sum of the indexes of the fixed points within it. 

Remark 2.3. Robustness of Poincare Index: 
The Poincare index is stable under certain deformations, 
the circle's Poincare index remains +1. 



4 



/ / \ \ 

(a) 



Figure 2: Robustness of Poincare index: start with any point on the circle and traverse circle counterclock- 
wisely. Upon returning to the starting position, the velocity vector rotates continuously for 2ir. Therefore, 
Poincare index of the circle is +1 in both cases. 




(a) (b) 



Figure 3: Two vector fields and their solution curves. In (a), the square is an isolating block. In (b) the 
square is not, because p is on the boundary and the solution curve passing through p is completely contained 
in the square. 

2.4 Conley Index 

Conley index is a powerful tool in dynamical systems. A fundamental theorem states that for a non- 
degenerate dynamical system, one can infer the structure of invariant sets within an so-called isolating 
blocks ./V by only studying the vector field along the boundary dN. Specifically, one computes the Conley 
index of N. Moreover, the index is invariant under certain deformation of the vector field. In image pattern 
recognition, one important task is to identify invariance of objects from a certain class. In many cases, these 
invariance are also invariant sets in related dynamical systems. However, due to various factors, the intensity 
landscape can deform or perturbate from sample to sample. Therefore we need a robust method to detect 
these invariant sets. We believe Conley index is a possible solution. Next, we review basic concepts and 
results in Conley index theory. For more details, see lfTTl l2l. 

Definition 2.4. A compact set N is an isolating block iff every boundary point leaves N immediately in 
forward time or backward time. 

See Figure [3]for example. 

Definition 2.5. Let X be the phase space with flow 4>. Let N be a compact set. The exist set N~ is defined 
to be 

N- := {x € dN : (f>(x, t) £ N Vt > 0} 

Definition 2.6. Let N be an isolating block, and N~ be the exit set. The Conley index is defined to be the 
homotopy type of the quotient space N/N~. 

Example 2.7. In planar flow, the Conley index for a source is S 2 , the 2-dimensional sphere. To see this, take 
N to be a disk around the source, then the boundary circle dNis the exit set N~. Recall that the quotient 
space N/dN is obtained by identifying dN to a point, which gives a sphere. Cf. Figure @ One can show 
that Conley index for a saddle is S , and a sink has index S°. For details, see £77]/. 
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Figure 4: Conley index for a source 



3 Hamiltonian Streamline Based Features 

Given an image I(x, y), we can derive a negative gradient flow system and a Hamiltonian system based on 
/. We believe that solution curves of the Hamiltonian system and invariant sets in the gradient flow system 
are important features of the image. For example, in face detection, eyes can be identified as attracting sets 
and nose as repelling sets. In this section, we propose four features, which are built based on orbits in the 
Hamiltonian system induced by image landscape I(x, y). Let o be an orbit in the Hamiltonian system, we 
have 



1. Density match. Let d(o, I) = I \ a , i.e., the density landscape restricted on orbit o. If o is of length n, 
then d{p, I) can be identified as a point in W 1 . For any other image /, one can measure the similarity 
between I and I by s(o; I, I) = | \d(o, I) — d(o, I)\\2- Suppose / and I are images from the same 
class, say human face, and o is a feature, say the boundary of an eye. We would expect s(o; /, /) to 
be relatively small. In contrast, if / and I are from different classes, we expect s(o; /, I) to be large. 
A limitation of density match is that it may not be robust under change of lightning. Suppose / and 
I are human face images under different lighting conditions, then their image landscapes may have 
different overall latitude, this may result in a large s(o; /, I) value. 

2. Direction match. Let 




be a negative gradient flow system. By normalizing this vector field, one can associate the image I 
with a direction field, vi(x,y). One can identify the direction field as a map <pj : M 2 — » M, where 
<p(x,y) = arctan(^). Let dir(o,vi) = 4>i \ . For an orbit of length n, dir(o,vj) G W 1 . One can 
define a similarity function, such that for any given image /, s(s; I, J) = | \dir(o, vj) — dir(p, f/)||2- 
Note, the negative gradient flow along an orbit describes relative latitude of nearby landscape around 
o. Under uniform lightening change, image landscape may shift up or down, but overall shape stays 
similar. In this case, direction match gives a more robust result than density match. 



Moreover, if o is a closed orbit, we have additional features 

3. Poincare index of o 

4. Conley index of the compact region enclosed by o 

In the next section, we will describe the methods to extract orbits in the Hamiltonian system and how to 
compute the Poincare index and Conley index numerically. 
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Figure 5: Forward algorithm 

4 Numerical Treatments 
4.1 Solution Curves 

The very first task in our method is to extract solution curves. However, there are two difficulties: (1) the 
dynamical system we obtained is usually a very sparse sample of an underlying "smooth" vector field, hence 
long time numerical integration may cause large error; (2) we must restrict solution curves on the lattice, 
therefore traditional numerical methods do not directly apply. 

Since solution curves are restricted on the lattice, the question is: given current position, which of the 
neighborhood points should be considered as the next point along the solution curve? Here, we propose a 
second order Runge-Kutta like algorithm. The idea is that we first integrate half step size to obtain a middle 
point; then determine the final position by estimating velocity field at the middle point. In particular, we 
estimate the velocity at the middle point as a weighted average of velocities at nearby points. Cf. Figure 
[5] This middle-point scheme is reminiscent of the second-order Runge-Kutta method, which has better 
accuracy than Euler's method. Empirically, our method gives satisfactory results. 

Algorithm 1 Forward Algorthim 

1: Input: vector field v in x-y coordinate and current position in image coordinate. 

2: Output: image coordinate of the next point along solution curve. 

3: mid-point z [j, i] + 0.5 x v(j, i) 

4: Let pi , p2 , P3 , Pi be the four nearest grid points. Set v, <— v \ Pi . 

5: Assign weights Wi oc ,, 1 i = 1,2,3,4 

\\Pi z \\ 

6: Set v z <- J2i Wi^i 

7: Let p* be the vertex such that the angle between p* — z and v z is minimized. Output the image coordinate 
of p*. 



The next question is how to generate a complete orbit based on our forward algorithm. Let 4>t De a flow 
and xq be an initial position. Let n be a non-negative integer. By (j) n (xo) we denote the point obtained by 
following the flow for n unit time. By 0_ n (xo) we denote the point obtained by following the backward 
flow for n unit time. A complete orbit on the lattice can be represented by {4>i(xo)}i, where — m < i < n 
for some non-negative integer m, n. The forward algorithm discussed previously provides an approximation 
for the map <p\. By applying the forward algorithm for the backward flow, we have an approximation for 
<fi-i. Note that 4> n = {4>i) n and = (</>_i) n , this gives us an iterative way to generate the complete orbit. 
However, we stop either we hit the boundary of the image or we detect a self-intersection. In Figure |6l we 
show three orbits generated by our algorithm for the Hamiltonian system in Figure Q] In Algorithm |2j we 
show a method to generate all solution curves for a system. 
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Algorithm 2 Extracting all solution curves 



1: Let S be the set of all known solution curves. Initialize 5 = 

2: Start with any initial condition that is not on any existing solution curve in S. 

3: extract the complete orbit o 

4: S=SUo 

5: Repeat 2-4 until we have exhausted all possible initial conditions. 




Figure 6: Some orbits in the Hamiltonian flow generated by face image. These orbits successfully capture 
boundaries of eyes, nose, forehead, cheeks, mouth, etc. One should note that the dataset is small, so the 
intensity landscape of the average face is not entirely symmetric. This makes these orbits look imperfect. 

4.2 Numerical Treatment of Poincare Index 

Recall that the Poincare index is defined as an integration along a closed curve. Here, we are concerned with 
closed orbits from the Hamiltonian system and evaluate the Poincare index of these curves in the negative 
gradient flow system. The heuristic is that closed orbits from the Hamiltonian system usually enclose inter- 
esting regions, like sinks and sources of the associated gradient system. Due to the robustness of the index, 
we expect this to be a robust feature. 

Unlike continuous flow, when working with images, we have to do integration numerically. Let o be a 
closed orbit of length N + 1, with o(l) = o(N + 1). Let v(x, y) be the normalized gradient system. Then 
the direction field along o is (f> = v(x, y) \ Q . Here, each 4>{i) represents the direction of the vector field, and 
we can identify them as angles in [0, 2ir). With this in mind, we compute the Poincare index as 

1 N+l 

_£(0(i)-0(i-l)) 

i=2 

Evaluation of the above equation is straightforward. However, there is a subtlety: let a 2 and a\ be two 
consecutive angles along the integration path, and we want to evaluate a 2 — a\ , then 

{a 2 + 2vr - a x if a 2 € [0, tt/2], ai G (3vr/2, 2vr) 
a 2 - (al + 2vr) if a 2 £ [3vr/2, 2tt), a% £ [0, vr/2) 
a 2 — a\ otherwise 

Another thing to note is we assume the closed orbit o to be counterclockwise. This can be done, for 
example, by built-in matlab function poly2ccw(). However, this is not absolutely necessary. If o is clock- 
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(a) D 2 (b) S° (c) S 1 (d) S 1 V S" 1 (e) S 2 

Figure 7: Examples of isolating blocks and exiting sets for planar flow. Red edge represents a exiting set. 
By identifying exiting set to one point, we obtain a quotient space whose homotopy type is the Conley 
index. (View in color). 

wise, then the index computed will be the actual Poincare index multiplied by —1. This will not affect the 
effectiveness of the feature, as long as one convention is followed. 

4.3 Numerical Treatment of Conley Index 

In general, it is highly complicated to compute the Conley index, because the index is the homotopy type 
of a quotient space. For those who are interested in the general method, we recommend Q. However, 
in this paper, we are only interested in planar flow, and we restrict our isolating blocks to be disks (up to 
homeomorphisms). One can show that under these assumptions, the Conley index is one of the following 
types: D 2 ,S°, S 2 and V^ 1 S 1 for some N > 1. Cf. Figure |7J See O for details. Instead of computing 
the index rigorously, we propose two alternative approximations, which we call continuous pseudo Conley 
index and discrete pseudo Conley index. 

1. continuous pseudo Conley index: Given an isolating block N and its boundary F = dN, Let N~ 
be the set of exiting pixels, i.e., pixels at which the vector field pointing outwards. We define the 
continuous pseudo Conley index to be card(N~)/card(F). Note, if the true index is 5°, this ratio is 
0; if the true index is S 2 , this ratio is 1. However, if the true index is S l , such pseudo index has 
no virtual connection with the true index. 

2. discrete pseudo Conley index: Let N be the number of connected exiting components on the boundary. 
Then the discrete pseudo Conley index is 

(i) D 2 , if there are exactly one connected exiting component and one connected entering component 

(ii) 5°, if there is no connected exiting component 

(iii) vf^S 1 if there are N > 2 connected exiting components 

(iv) S 2 if there is one connected exiting component and no pixel is an entering pixel 

In this work, we choose to compute the continuous pseudo Conley index, because this pseudo index 
recovers the true index in case of attractor and repellor. We are more interested in regions of attracting and 
repelling, and these index serve our purpose sufficiently. 

5 General Framework 
5.1 Canonical image 

Our method relies on the hope that we can actively extract possible features of a given object class. Such 
features are built on two dynamical systems induced by the image landscape. However, it is not favorable 
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to work on a single sample image, since features extracted this way may not be representative. Instead, we 
prefer to work on a canonical image, which provides a "standard" image landscape for this object class. We 
extract feature lines from this canonical image, and use them as feature templates. We expect that images 
from the same class have similar feature values. 

One question is how to obtain such a canonical image. In general, this can be highly nontrivial. However, 
if the object class has non-deformable shape and training images are well aligned, one can use the average 
image as a canonical image. In Figure [H we show the average image and a random sample in our training 
set. The landscape of the average image is more desirable than that of the random sample. 

Our general framework has two phases: I) Design a family of features based on a canonical image. II) 
Use these features as weak classifiers and apply AdaBoost ifTTl to obtain a final strong classifier. We refer 
readers to IfTTl for the AdaBoost algorithm. Here, we give a brief outline for phase I): 

1. Start with a canonical image and compute its negative gradient system and Hamiltonian system. 

2. Extract all orbits in the Hamiltonian flow. 

3. Any orbit is associated with at least two features: density feature and gradient direction feature. 

4. If an orbit is closed, we associate two more features with it: the Poincare index and the Conley index. 

5. We evaluate feature values according to their types. 

5.2 Experiment Results 
5.2.1 The Datasets 

We have done experiments on the AT&T ORL dataset and the MIT CBCL dataset. (Thanks to AT&T 
Laboratories Cambridge and MIT CBCL). The ORL dataset contains 400 face images of 40 people. Images 
are of size 112 x 92. To get negative samples, we randomly extracted 800 non-face images from the Caltech- 
256 O clutter set. We use 200 faces and 400 non-faces for training and the rest for testing. The MIT CBCL 
dataset contains images of size 19 x 19. Training set contains 2429 faces and 4548 non-faces, the test set 
containing 472 faces and 23,573 non-faces. 

For both datasets, we compare our Hamiltonian classifiers with Haar-like classifiers. A comparison with 
SVM is also included. 



10 





Hamiltonian 


Haar 


10-stage classifier 


1 fn, fp 


3 fn, 4 fp 


20-stage classifier 


fn, fp 


2 fn, 1 fp 


100-stage classifier 


not needed 


1 fn, fp 



Figure 9: fp = false positive, fn = false negative 



5.2.2 Experiment on the ORL Dataset 

We take the average of training faces, and use it as a canonical image. The intensity landscape of the model 
image gives us 45 orbits in the associated Hamiltonian flow, of which 20 are closed. This gives us 130 
Hamiltonian streamline based features. We also generate Haar-like features of those types suggested in 
ITP71 . For training images of this size, there are about ( 9 2 2 ) ( 1 2 2 ) 2.6 x 10 7 possible Haar-like features. 
However, in our experiments, we uniformly choose about 27,000 such features. 

Using Hamiltonian Streamline Features If we only use the 130 Hamiltonian streamline features, training 
time is very short: approximately 2.5 minutes on our machine (Intel(R) Core(TM)2 Duo CPU P8600 @ 
2.40GHz 2.40GHz). Figure [6] shows the first few features selected by AdaBoost. 

Since we have designed four types of Hamiltonian streamline features, one natural question is in what 
order these types of features are selected in AdaBoost. In our test, we find "direction match" is always 
the first selected feature, followed by "density match". Conley index comes third, while Poincare index is 
least favorable. In a 20-stage classifier, 6 of which are "density match", 1 1 of which are gradient "direction 
match", 1 of which is Poincare index, 2 of which are Conley index. 

Comparing with Haar-like Features One would wonder how our Hamiltonian streamline features stand 
against Haar-like features. If we use 27,000 Haar-like features, AdaBoost takes more than 1 hour on the 
same machine. Testing accuracy is slightly worse than the Hamiltonian classifier. See Figure [9] 

Joint selection: Hamiltonian streamline features and Haar-like features It would be interesting to 
see, when pooled together, which type of feature is favored by AdaBoost algorithm: Hamiltonian streamline 
features or Haar-like features? We put 130 Hamiltonian features with 27000 Haar-like features and let 
AdaBoost do feature selection. In Figure [10] we show the first 20 features selected by AdaBoost and their 
weights in the strong classifier. It turns out that Hamiltonian streamline features are more favorable. 

Detection Speed We have seen that with our Hamiltonian features, training speed is extremely fast. It is 
natural to ask whether the Hamiltonian features can achieve high speed in detection. The computational cost 
of evaluating Hamiltonian features is proportional (with a low constant) to the length of the orbits. In our 
experiment on the ORL dataset, the mean length of the orbits is 149 (median is 77). The 10-stage classifier 
has about the same accuracy as the 100-stage Haar-like feature classifier. In the 10-stage Hamiltonian 
classifer, total length of orbits is 1774. To evaluate a single Haar-like feature, we need to access 8-16 
pixels in the integral image (depending on its type). In total, we need to access 800-1600 pixels. The total 
operations needed in two methods do not have a significant difference asymptotically. Empirically, we find 
both classifiers have about the same speed in classifying one subwindow. Since Haar-like classifiers can 
do face detection in real-time, our feature can do as well. Moreover, boosting classifiers can be trivially 
parallelized; with today's GPGPU, real-time detection is not an issue. 



11 



Hamiltonian 
Haar 



o o o 



Figure 10: First 20 selected features, their types and weights in the final strong classifier 



5.2.3 Experiment on the MIT CBCL Dataset 

Since the ORL dataset is relatively small, we have also done experiments on the MIT CBCL dataset. Just 
like the ORL dataset, we have done 3 groups of experiments: Hamiltonian features alone, Haar-like features 
alone, Hamiltonian features plus Haar-like features. ROC curves are shown in Figure [TT] For this dataset, 
Haar-like classifier performs slightly better than Hamiltonian classifier. We believe this is due to the low 
resolution of training images whose underlying vector fields are not smooth enough. If we combine two 
types of features, the hybrid classifier achieves best accuracy. The ROC curve of the hybrid classifier is 
almost the same as the ROC curve of a SVM classifier using 2nd degree polynomial kernels 0. 

5.2.4 Discussions 

We summarize some observations in our experiments and discuss several advantages of Hamiltonian features 
over Haar-like features as follows: 

(1) Training with Hamiltonian features is much faster. Computational cost of an AdaBoost algorithm is 
proportional to the number of weak classifiers (features). For n x n training images, there are at least 
0(n 4 ) possible Haar-like features. However, empirically, our feature is of order 0(n). 

(2) Hamiltonian features are more descriptive and concise. Unlike Haar-like features which are ran- 
domly generated, Hamiltonian streamline features can intelligently capture class-related features. In 
an AdaBoost classifier, one needs far fewer Hamiltonian weak classifiers than Haar-like feature based 
classifiers. This helps to keep detection speed comparable with Haar-like classifiers. 

(3) For training images with higher resolution, Hamiltonian features have better performance. For low 
resolution training images, Hamiltonian features may perform slightly worse, but can be used jointly 
with Haar-like features to achieve better accuracy than using Haar-like features alone. Since the 
Hamiltonian features are built based on two vector fields derived from an image, their performance 
depends on how well we can approximate these underlying vector fields using algorithms in Section 
4. Algorithm 1 is a second order Runge-Kutta like algorithm, and the second oder Runge-Kutta 
method has error bound 0{l/h 2 ), where h is the step size. The ORL images have about 5 times more 
resolution than CBCL images (in one dimension). Thus, the step size for ORL is about 1/5 of that for 
CBCL. According to the error bound, orbits generated from ORL images are much more accurate. We 

1 SVM result obtained from |http ://vision.ai.uiuc. edu/mhyang/f ace -detect ion- survey ■ html| 
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Figure 11: ROC curves on the MIT CBCL dataset 

conclude that Hamiltonian features can achieve higher detection accuracy, if higher resolution images 
are used. 

6 Limitations and Future Work 

In this section, we discuss some limitations of our method and possible future works. Our method needs 
a canonical image. In our experiment, we take the average of training images as the canonical image. 
However, if training samples are not well aligned, the resulting average may be undesirable. In addition, if 
the object class has various possible shapes (e.g. cars), then it is not reasonable to use only one canonical 
image. One possible solution is to sample sufficient number of training images and use them as a family of 
canonical images. Then extract all features from all these images. 

In addition, our method is not robust against moderate rotations, because the first two feature are not so. 
However, since our training time is short, one can rotate the canonical image for a series of degrees and train 
a series of classifiers. But this will certainly increase detection time. 

7 Conclusion 

This paper presents a feature extraction method based on dynamical systems induced by image landscapes. 
Such features capture both topological and geometrical characteristics of the image landscape of an object 
class, and can be used in a boosting algorithm for object detection. In our experiments, we find training 
with Hamiltonian features is much faster than training with Haar-like features. In addition, for training 
images with good resolution, Hamiltonian features also yield better accuracy. In a final strong classifier, one 
needs much fewer Hamiltonian weak classifiers than Haar-like weak classifiers. We find detection speed of 
a Hamiltonian classifier is comparable to that of a Haar-like classifier. We believe our feature can be used as 
a good complement for existing features in the computer vision community. 
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